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Abstract 

This paper deals with the kinetic theory modeling of crowd dy- 
namics with the aim of showing how the dynamics at the micro-scale 
is transferred to the dynamics of collective behaviors. The derivation 
of a new model is followed by a qualitative analysis of the initial value 
problem. Existence of solutions is proved for arbitrary large times, 
while simulations are developed by computational schemes based on 
splitting methods, where the transport equations treated by finite dif- 
ference methods for hyperbolic equations. Some preliminary reason- 
ings toward the modeling of panic conditions arc proposed. 

Key words: Crowd dynamics, complexity, scaling, living systems, non- 
linear interactions. 

1 Introduction 

The modeling of pedestrian crowd dynamics can be developed at different 
representation scales, namely micro- and macro-scales with the intermediate 
approach of the kinetic theory, where the dependent variable is a probabil- 
ity distribution over the micro-state of pedestrians. The interested reader is 
addressed to the review paper [6], where the existing literature is reviewed 
and critically analyzed focusing at different modeling scales. Namely, at the 
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micro-scale by ordinary differential equations, see among others [18, 20, 21]; 
at the macro-scale obtained by the classical approach by conservation equa- 
tions, among others [10, 24]; or by stochastic models related to evolving prob- 
ability measure [26]. The interplay between individual-based and macro- 
scopic models is studied in [13]. Additional bibliography and conceptual 
links with the modeling of vehicular traffic are proposed in [5]. A detailed 
analysis of empirical data is delivered by various papers [9, 27, 30] among 
others. 

This present paper specifically refers to [3], where the main hint, fur- 
ther stressed in [6], is that the modeling approach should retain, as far as 
it is possible, the complexity features of crowds to be viewed as a living, 
hence complex, system. This paper presents an introduction to modeling 
where pedestrians move with only velocity modulus and change their ve- 
locity directions both due to interactions with other pedestrians and to the 
search of their specific target, for instance the exit zone. This present pa- 
per shows how the dynamics at the micro-scale of pedestrians is transferred 
to the dynamics of collective behaviors. Modeling and simulations include 
interactions between pedestrian that move toward different targets, while 
their dynamics can be induced, by signals from the outer environment to an 
exit direction. A qualitative analysis of the initial value problem, and some 
sample simulations are presented. 

The approach developed in this paper is based on the methods of the 
kinetic theory of active particles, which shows the ability to retain various 
complexity features. This approach has been applied to the modeling of 
several living systems as documented, among others, in the modeling of 
vehicular traffic [4, 8], swarms [7], social competition [1], and migration 
phenomena [23]. 

The contents are proposed through four additional sections. In details, 
Section 2 presents the mathematical model obtained as a development of 
that given in [3] to include the dynamics on the velocity modulus and inter- 
actions of pedestrians that move toward different directions [2]. Section 3 
introduces the problem of existence of solutions to the initial value problem 
in unbounded domains and announces the existence theorems for arbitrary 
large times. The technical proofs are developed in Sections 4 and 5. Sim- 
ulations obtained by application of the splitting method, are presented in 
Section 6 with a detailed computational analysis of the dynamics of the 
moving boundary of the domain occupied by the crowd. Finally, the last 
section first critically analyzes the use of the kinetic theory approach and 
subsequently looks ahead to research perspectives focusing on the modeling 
panic conditions. 
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2 A Mathematical Model 



Let us consider a large system of pedestrians moving in a bounded domain 
Q C R 2 , which may include, as shown in Figure 1, internal obstacles. The 
boundary of f2 is denoted by d£l. The modeling approach also refers to 
unbounded domains in two space dimensions. Pedestrians can be subdivided 
into two different groups, moving toward two different targets, denoted by 
Ti and T2, respectively. 




Figure 1: - Geometry of the domain fi: Pedestrians enter from I and can 
move towards two different targets: pedestrians at P and at Q move to T2 
and Ti, respectively. 

A preliminary step of the derivation of the model, following [3], consists 
in the assessment of the representation of the system according to the meth- 
ods of the generalized kinetic theory. This means selecting the microscopic 
state of the individual state of the pedestrians and of the probability dis- 
tributions over such state. Two distributions are considered, corresponding 
to the aforesaid groups of pedestrians pursuing different objectives, which 
are the dependent variables of the equation modeling the overall dynamics 
of the system. 

Bearing all above in mind, the microscopic state can be defined by the 
variables: position x = {x,y} e!]C i? 2 , and speed v 6 D v C R 2 . Dimen- 
sionless coordinates are used by dividing space coordinates by the largest 
dimension I of O, while the velocity modulus v is referred to the largest ad- 
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(a) (b) 

Figure 2: (a) Interaction domain for an individual in the position x. (b) 
Interaction between a candidate particle in the position P moving in the 
direction 9^ and a field particle in Q moving in the direction 9 P ; 9 V is the 
angle from P to the target T. 

missible velocity vi that can be reached by a speedy pedestrian in optimal 
conditions of the ambient where the crowd moves. 

Polar coordinates with discrete values are used for the velocity variable 
v = {v,9}, where v is the velocity modulus and 9 identifies the velocity 
directions with respect to the x-axis: 

Tl 

I e = {9i =0,... ,9i,... ,9 n = -2tt}, 

n — 1 

I v = {vi = Q,... , Vj , . . . , v m = 1}. 

Out of this discrete velocity representation some specific cases can be 
considered. For instance the case of two different groups, labeled with the 
superscript a = 1,2, which move towards two different targets. Moreover, 
if the activity variable is assumed to be homogenously distributed for all 
pedestrians, the corresponding representation is as follows: 

n m 

/ ff (t, x, v, u) = J2 E ffj (* , *) 8(0 - 9i) ® 5(v - Vj )® 5(u - u ) , (1) 

i=l j=l 
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where f[j(t, x) = f a (t, x, 9i, Vj, no) corresponding, for each group a = 1,2, 
to the r/-particle, namely to the pedestrian moving in the direction 0j with 
velocity Vj. 

Macroscopic quantities are obtained by weighted sums. In particular, 
the number density and flow of pedestrians are, respectively, given by: 

2 2 n m 

x) = Y /P a (t,x) = Y,Y,Y, fij (*. x ) . ( 2 ) 

(7=1 <T=1 1=1 j = l 

and 

2 2 n m 

q(t,x) = ^q ff (t,x) = cose i)Wj sin^)/5(i,x). (3) 

<T=1 (7=1 1=1 j = l 

Remark 2.1. The use of discrete space microscopic states was introduced 
in [15] and [12] to model the granular behavior of vehicular traffic, namely it 
is observed in [14] that the number of vehicles is not large enough to justify 
the assumption of a continuous distribution function. The same reasonings 
are valid also in the case of pedestrians in a crowd. 

Remark 2.2. The more general case, where the activity variable is hetero- 
geneously distributed and is modified by interactions is critically analyzed 
in the last section focusing on the modeling of panic conditions. 

The derivation of models needs to be related to a mathematical structure 
suitable to retain, as far as possible, the complexity of crowds dynamics. 
Accordingly, the following qualitative features of the dynamics are taken 
into account: 

i) Pedestrians have a visibility zone A = A(x), see Figure 2(a), which does 
not coincide with the whole domain f2 due to the limited visibility angle of 
each individual. The specific geometry of f2 can reduce the visibility zone 
for pedestrians near walls. 

ii) Interactions, at each at time t, involve: the test (ij')-particle, which is 
representative of the whole system; field (pg)-particles in x* of the visibility 
zone A; and candidate pedestrian (hk), in x. The candidate pedestrian 
modifies, in probability, its state into that of the test pedestrian, due to 
interactions with the field pedestrians, while the test pedestrian losses its 
state due to interactions with the field pedestrians. The activity of the 
candidate pedestrian is u w , which is not modified by the interaction. 

iii) The output of interactions is modeled, by using tools of the game the- 
ory [28, 29], by the discrete probability density function AJ lk pq (ij), which 
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denotes the probability that a candidate (/i/c)-particle modifies its state into 
that of the test (ij')-particle due to interaction with the field (pg)-particles. 
The transition density depends on the whole local number density p. 

iv) The frequency of the interactions is modeled by the term 77 that is allowed 
to depend on the density of the field particles in the visibility zone. 

Time and space dynamics of the distribution function can be obtained 
by equating its increase on time in the elementary volume of the space of 
the microscopic states to the net flow into such volume due to transport and 
interactions: 

(d t + Vtf ■5J/5(t,x) = ,7[f](*,x) (4) 
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where f = {fij}, while the term •Af lkpq (ij) should be consistent with the 
probability density property: 

n m 

EE^«(y) = 1 > V/ip€{l,...,n}, Vfcg€{l,...,m}, (6) 
i=ij=i 

for a = 1,2, and for all conditioning local density. The derivation of spe- 
cific models can take advantage of the above structures and is obtained by 
modeling the various terms that appear in Eq. (5). The modeling necessar- 
ily needs heuristic assumptions towards the interpretation of the complex 
phenomenology of the system under consideration. 

• Interaction rate: The modeling of the term 77, can be developed sim- 
ilarly to the case of vehicular traffic [8], by increasing the interaction rate 
with increasing local density in the free and congested regimes. For higher 
densities, when pedestrians are obliged to stop, one may assume that 77 keeps 
a constant value, or may decay for lack of interest. The following model can 
be proposed among various conceivable ones: 

V (p(t, x)) = 7?°(1 + p(t, x)) exp ( - p(t, x)) , (7) 

while a simpler model can be obtained by assuming 77 = 77°. 

• Transition probability density: The modeling of the transition prob- 
ability density -A% kpq (ij) refers to the interactions involving candidate and 
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field particles. The approach proposed here is based on the assumption that 
particles are subject to three different influences, namely the trend to the 
exit point, the influence of the stream induced by the other pedestrians, and 
the selection of the path with minimal density gradient. A simplified in- 
terpretation of the phenomenological behavior is obtained by assuming the 
factorization of the two probability densities modeling the modifications of 
the velocity direction and modulus: 

A a hk;Pq (ij) = Bl p {i){6 h -> 0i|p(*,x)) x Cl q (j){v k -> t)>(t,x)), (8) 

where subscripts correspond to the state of the interacting pairs, superscripts 
to the group and the output of the interaction is indicated in the argument, 
which depends on the local density. 

Remark 2.3. Equation (8) is based on an heuristic assumption that can- 
not be justified by theoretical issues. However, it simplifies the modeling 
approach by separating the various causes that modify the dynamics. This 
assumption does not imply factorization of the probability density due to 
the mixing action of the variable p. 

Bearing this remark in mind, let us consider, separately the modeling of 
the two probability densities: 

- Modeling B^ p (i) : Let us consider the candidate particle in the position 
P moving in the direction Oh and interacting with a field particle with di- 
rection p . Moreover, let V be the angle from P to T a or to the direction 
of the shortest path, see Figure 2(b). Let us introduce a parameter a a , 
which models the sensitivity of pedestrians both to reach the target T a , and 
to follow the stream. Moreover, it is assumed that the sensitivity to the 
target, depends on the ratio of vacuum 1 — p, while the sensitivity to the 
stream depends on the density p. Finally it is assumed that the dynamics 
is more active for higher values of the activity variable uq. It is natural 
that an active individual shows a higher ability to modify the direction of 
the motion. The games modeling interactions are the same for both groups 
being understood that the target chased by pedestrians differs for each of 
them. This qualitative description can be formalized as follows: 

- Interaction with a upper stream and target directions, namely P > 
9h, U > Oh: It is assumed that, in addition to the selection of the path 
with minimal gradients, both actions contribute to an anticlockwise rotation: 

Bh p (i) = a a u (l - p) + a a u p if i = h + 1 , 
Bhp{i) = 1 — Oo- uq(1 — p) — a a no p if i = h, 
Bl p {i) = if i = h-l. 
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Analogous calculations can be done for the other cases. The result is as 
follows: 

- Interaction with a upper stream and low target direction 9 P > 0^; 6 U < 

e h : 

Bl p (i) = a a u p if i = h + l, 

Bhpii) = I - a a u (l - p) - a a u p if i = h, 

®h P {i) = a a u o (1 - P) if i = h—l. 

- Interaction with a lower stream and upper target direction 9 p < 9^; 9 V > 

e h -. 

Bh P {i) = "o-^o(l-p) if i = h + l, 

®hpS) = 1 ~~ a <y u o(l — P) — ct a no p if i = h, 

®hpS) = ct a uo p if i = h — 1 . 

- Interaction with a lower stream and target directions 9 P < 9^; 9 V < 9^: 

Bl p {i) = if i = h + l, 

B% p (i) = 1 - a CT u (l - p) - a a u p if i = h, 

Bh P (i) = ct a u (l - p) + a a u p if i = h-l. 

Concerning the term C% q (j), the proposed modeling assumes that pedes- 
trians have a trend to adjust each other their velocity modulus, namely an 
encounter with a faster (slower) pedestrian induces a trend to increase (de- 
crease) the speed. This trend is conditioned by the local density, namely 
it decreases with increasing density. It is assumed that also this type of 
dynamics depends on the activity no of the pedestrian and on their sensitiv- 
ity (3 a to the change of velocity. The corresponding table of games, which 
defines the interaction rules, is as follows: 

• Interaction with faster (or slower) particles Vk < v q , or > v q , respec- 
tively 

!1 - P a u p, j = k; 
PaU p, j = k + l; 

0, otherwise; 
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Cl q {3) 



Pa U p, 
1- Pa U p, 
0, 



j = k; 
3 = k - 1; 
otherwise. 



• Interaction with equally fast particles k = q 

W = 



1 - 2/3 CT u p, 

Pa U p, 
Pa U p, 



j = k; 
j = k - 1; 
j = k + 1. 



Moreover, 

• for k = 1 i/ie candidate particle cannot reduce the velocity, while for k = k 
cannot increase it: 



c a kq {j) 



Cl q {j) 



1 - P a U p, 
Pa U p, 

o, 

Pa U p, 

1 - Pa U p, 

o, 



i = i; 

3 = 2; 
otherwise; 



J 



m 



1; 



j = m; 
otherwise. 



• The modeling corresponding to jam densities p G [e, 1] implements that 
pedestrians reduce to zero their velocity. 

Remark 2.4. If the heterogeneous behavior of pedestrians needs to be 
taken into account a probability distribution over the variable u has to be 
introduced with the modeling of the interaction dynamics that modifies it. 



3 On the Initial Value Problem 

The initial value problem consists in solving Eqs. (5) with initial conditions 
given by 

/5(0,x) = ^(x). (9) 

Let us introduce the mild form obtained by integrating along the char- 
acteristics: 

}g(t,x) = ^.(x) + ^(fg.[f,f]( S ,x)-^.( S ,x)L[f]( S ,x))d S , 

i€{l,...,n}, j e{l,...,m}, a €{1,2}, (10) 
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where the following notation has been used for any given vector /(i,x): 
f[j{t, x) = fij(t, x + Vj cos(0j)i, y + Vj sin(^)t), and where the operators in 
(10) are defined as follows: 



n m ,. 

ry.[f,f] (t,x)=£ 2 / ^P(i,x*)]^(u)[p(t,x*)] 

h,p=l k,q=l A 

/ftfc (*> x + ( v j cos(0i) - v k cos(6 h ))t, y + (vj sin(0j) - w fc sin(^))t) 
/pq - v q cos(e p )t,y* - v q sm(6 p )t) dx* , (11) 



and 



L[f](t, x) = ^£ f7[p(t,x*)]/^(t,x*- V ,cos(0 p )t,y*-t;,sin(e p )t)dx*, 

(12) 

where f = {/^}- 



Let us now define the following functional space: 

2 n m 



rl 



{/ = (/5)eM 2n>2m :||/|| 1 =^^^ / |/5-(t,x)| dx<oc} 

CT =1 j = l j = l Jn 

and introduce, for T > 0, the Banach space 

X T = X[0,T] = C([0,T],L^ 2m ) 

of the matrix-valued functions / = /(i,x) : [0, T] x 17 — >■ Mm^m- such that 
for all fixed t £ [0, T] the function x — >■ /(i,x) belongs to -^M 2n2 m' wnere 
the norm in Xt is given by 

|| / \\x T = sup || / 111 . (13) 

te[o,T] 

The qualitative analysis developed in the following takes advantage of 
the following assumptions that are consistent with the physics of the system 
under consideration: 

Assumption H.l. For all positive R, there exists a constant C v > so 
that 

< rj(p) < C v , whenever < p < R . (14) 

Assumption H.2. Both the encounter rate rj[p] and the transition prob- 
ability A a hk pq(ij)[p] are Lipschitz continuous functions of the macroscopic 
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density p, i.e., that there exist constants L v ,La such that 

I V[Pl] -V[P2] \< L v \p!-p 2 |, 

(15) 

I ^fe,pg(u)[pl] - ^fe,pg(u)[P2] |< ^ | Pi ~ P2 |, 

whenever < p\ < R, < p2 < -R, and all i,h,p = l,..,n and j,k,q = 
1, ..,m. 

The following results can be proved under assumptions H..1 and H.2: 

Theorem 3.1 (Tocal Existence,) Let cp^ £ Loo n L 1 , ^ > 0, £/ien £/iere 
exists (f)° so that, if \\ <f> ||i< 0°, i/iere exist T, a®, and R so that a unique 
non-negative solution to the initial value problem of Eqs.(5)-(9) exists and 
satisfies: 

f G X T , sup || f{t) ||x< a || ||x, (16) 
te[o,T] 

and 

p{t,-x.)<R, Vt€[0,T], xeft. (17) 

Moreover, if 

2 71 171 

EEEii ^iu<i, as) 

cr=l i=l j=l 

and || ||i is small, one has 

p(t,x)<l, Vt€[0,T\, xefl (19) 



Theorem 3.2 (Global Existence,) Under the assumptions of Theorem 
3.1, there exist (j) r , (r = l,...,p — 1) such that if \\ 4> ||i< (ff , there ex- 
ists a r so that it is possible to find a unique non-negative solution to the 
initial value problem of Eqs.(5)-(9) satisfying for any r < p — 1 

f(t)eX[0,(p-l)T], (20) 

sup || f(t + (r-l)T) ||i<a r _i || ||i, (21) 
te[o,T] 

and 

p(t + (r-l)T,x)<12, Vte[0,T], xe!l. (22) 
Moreover, if (p satisfies (18), one has 

p(t+ (r- l)T,x) < 1, VtG[0,T], xe!1. (23) 
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4 Proof of Theorem 3.1 (Local Existence) 

Let us introduce the function 

^(t,x) = exp(At)/5(t,x), for A>0. (24) 

Therefore, system (5)- (9) can be written in equivalent form in terms of the 
functions ip(t,x) = (ip£j(t,x)), where ipfj are solutions of the following mild 
problem: 

^■(t,x) = ^(x) + | o *^K^(^x)exp(-A S ) 

+ ^.( a ,x)|A-L^]( S ,x)exp(-A S )})da. (25) 
Moreover, let us consider the operator A acting on Xt whose components 



are 



AtyO£(t,x) = ^.(x) + ^(rf J .[^V](5,x)exp(-A S ) 

+ ^.( S ,x){A-L^](s,x)exp(-A a )})cfo. (26) 

The principle of contractive mapping can be applied to Eq.(26). How- 
ever, the proof is based on various sharp inequalities that are preliminary- 
proved in the following Lemmas. 

Lemma 4.1 Let T > 0, ip G X T , then Atp £ X T and 3C\ > so that 

II A(ip) \\ Xt <\\ <t> Hi +y II $ tx T +AT || ^ \\x T ■ (27) 

Proof of Lemma 4.1: To prove the estimate (27), we need first to estimate 
the nonlinear operators Ffj[ip,ip] and L[ijj] in L\. Taking into account the 
assumption H.l on 77 yields 



r*#,v]||i < c v ]T E [ Ahkw&Mt'**] 

h,p=l k,q=l JA 

x | -tjjg q (t,x* — v q cos(6p)t,y* — v q sm(9 p )t) | dx* 

x / I + (wjcos(^) - v k cos(6 h ))t,y 

Jn 

+ (vjsm(0i)-v fc siii(0/,))*) | dx. (28) 
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Using Eq. (6) on the probability density A^ kpq (ij) and the following 
transformations with jacobian equal one 

x = (x, y) i-)- (x - v q cos(0 p )i, y - v q sin(0 p )i), 

x = (x, y) (->■ (x + (uj 003(6*1) - v fc cos(0h))t, y + (v,- sin(0j) - v k sm(6 h ))t), 
yields: 

II II 1— E E II rg#,Y>] ||i< G, || ^ II? • (29) 

cr=l i,j 

Moreover, the same estimate can be obtained for the loss operator L: 

II Lfi>$ lli< ^7 II V> II? • (30) 

The following inequality is directly deduced, by using (29)-(30), from 
(26): 

II Aty) ||i<|| ||i +^(1 - exp(-At)) || $ \\ 2 Xt +Xt || V ||x T , (31) 

which gives (27) with Ci = 2C^. 

Analogous calculations lead to the following inequality: 

|| Aty)(t) - A^)(s) ||x<| t - a I (d || V ||x T +A || V> ||x T ) , 

which gives the continuity of the operator A in L\. Therefore, Lemma 4.1 
is proved. ■ 

Lemma 4.2 Let i> 1 = (ipfj 1 ) G Li, V 2 = (V>f/) G £1, such that 

2 n m 

EEE^ 1 ^' x - vj cos(9i)t,y - vj sm(6i)t) < Rexp(Xt), 

a=l i=l j = l 

(32) 

2 n m 

2 S S^ 2 (*' x ~ u i cos (^)*' 1/ - u i sin (^)*) < Rexp(Xt), 

{ a=li=lj=l 

forxen,t> 0. Then 

|| f^Y^-r^ 2 ^ 2 ) 111 < (C v + nmR(C ri L A + L r! )) 

x (|| V 1 111 + H 2 111) II V' 1 -^ 2 111, (33) 

IIVAA^-V^A^ 2 ) ||i< {C v + L n R){\\ y, 1 ||i + || V 2 111) || V' 1 -^ 2 Hi • 

(34) 
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Proof of Lemma 4.2: First write T^ip 1 , ip 1 )^, x) - fg^ 2 , V> 2 )0, x) 
as follows: 



ram. 

= E E / ^M*' x *)]^-W^)[M*> x *)l 

U 1 7. 1 J A 



h,p=l k,q=l ' 



x V ; L 1 (*' X + (ujcos(^) - v k cos(0 h ))t,y + (uj sin(0j) - u fc sin(^))t) 
xil)°q{t,x* - v q cos(6 p )t,y* - v q sm(0 p )t) 

-^hk^i x + ( u i cos (^) ~ v k cos(9 h ))t, y + (vj sin(6»i) - t> fc sin(^))t) 



x^ q 2 (t, x* - v q cos(6 p )t, y* - v g sin(9 p )t) 



pq 

n m 



+ E E / K,k 2 (^ x + (vjcos(0i) -v k cos(e h ))t,y 

h,p=l k,q=l A 

+(vj sin(6»i) - v k sm(e h ))t) 

x ^p, g 2 (*) x * ~ u g cos(6 p )t, y* - v q sm(0 p )t) 



r)[p^i (t, x*)]Al km {ij) [p^i (t, x*)} 

-v[p^{t,x*)]Ahk, pq (ij)[p^{t, x*)] dx* 
= A ij + B ij, 



Using Assumptions H.l, H.2, and summing over yields 

2 n m 

EEEH 

a=i i=i j=i 

2 n m n m ,. 

< c v E E E E E L „ I ^^a, x + («,- 

CT =1 i=l j=l h,p=l fc,g=l JAxf7 

-u fe cos(6 h ))t, y + (vj sin(0;) - w fe sin(0 h ))t) | 
x | ip^it, x* - v q cos(9 p )t, y* - Vq sm(9 p )t) 
—tp pq 2 (t, x* — v q cos(9 p )t, y* — v q sm(9 p )t) \ dx*dx 

2 n m n m „ 

+^EEE E E/ A a ^m 

a=l i=l j=l h,p=l k,q=l JAXU 

x | i/>° 2 (t, x* - Vq cos(9 p )t, y* - v q sm(9 p )t) I 
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x I V^fc 1 {t, x + (vj cos(6»i) - v k cos(9 h ))t, y 
+(vj sin(^) - w fc sin(0 h ))i) 

~^h,k 2 (^ x + ( w i cos (^) ~ u fe cos(9 h ))t, y 
+(vj sin(0j) — Ufe sin(0/ l ))t) | dx*dx 

< C„ || ^ - V 2 111 (|| Y* 1 ||i + || V 2 Hi) ■ (35) 

On the other hand, considering that 

2 n m 

I (PV 1 ~ PV 2 )(*' X *) I- ex P(~^) V>ri/(*> ^* - V y COs{e r )t,y* 

c=\ r=l u=l 

—v u $m{Q r )t)iJ)° u 2 {t, x* — v u cos(6 r )t,y* — v„sm(9 r )t) \ , 

and 

n m 

E E I %?tt> x * " v q cos(O p )t,y* - v q sm(6 p )t) |< i?exp(At), (36) 

p=l g=l 

then, using Assumptions H.1-H.2, noting that ^ kpq (ij) < 1, and sum- 
ming over i,j,a, yields: 

2 n m 

E E E H B ij II + C„L A ) || V 2 || ill Y* 1 - V' 2 111 

(7=1 1=1 J = l 

which gives, by using (35), the estimate (33). 

By the same arguments, from Assumptions H.l, and H.2 one deduces 
(34) and this completes the proof of Lemma 4.2. ■ 

Lemma 4.3 Let T > 0, i> x = (Vf/),^ 2 = (Vf/) G Ay, satisfying (32) for 
x € fi,t € [0,T]. T/ien: 

L There exists a constant C2 > suc/i i/tai 

M(^)-A(^)||x T < (^(11^11^ + 11^11^) + AT) 

II V' 1 - Y> 2 ||x T • (37) 
2. If (V" 1 )^' > 0; ^ere 3Ao such that if A > Ao and i/<% > one has 

(A&m > 0. 
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3. Let (f)fj G L°° , then for R large enough(R > R\), there exists T such that 

2 n m 

E E E W 1 ))^, x - Vj cos(9i)t, y - Vj sm(6i)t) < Rexp(Xt), (38) 

(7 = 1 1=1 j = l 

Vt€[0,T\, V^gO. 

Proof of Lemma 4.3: From (26) and Lemma 4.1, one gets easily (37) with 
C 2 = C v + nmR(C ri L A + L v ). (39) 

Since ff^ip 1 , ip 1 ]^, x) > if (Y* 1 )^- > 0, the nonnegativity of (A^))^ 
depend on the possibility to find A > such that 

A - Lf^Ks, x) exp(-As) > 0. (40) 

Noting that from (12) and assumption H.l, one obtains 

n m „ 

L[^](t, *) = EE / r)[p(t,x*)]^; q (t,x* -v q cos(e p )t,y* 

p=lq=l JA 

-v q sm(9 p )t)dx* < Rexp(Xt)C v | A | . 

The nonnegativity of (A^ 1 ))^ is then achieved by choosing A > Ao = 
RC V | A | . To deal with 3) first let see that 

r(^ 1 , il> 1 )(t, x — Vj cos(9i)t, y — Vj sin(^j)t) 

= E E / r,[p(t,x*)]AZ k!Pq (ij)[p(t,x*)} 

h,p=l k,q=l A 

x^l ^t, x-v k cos(9 h )t, y-v k sm(6 h )t) 
x-0 1 p i g(i, x* — v q cos(9 p )t, y* — v q sm(6 p )t)dx* 

< C V R 2 exp(2At) f dx* , < C v | A | R 2 exp(2Ai). (41) 
J A 

Then the following estimate is proved: 

2 n m 

E E E (^W0)&(*> x ~ v i cos (^ V - vj am(0i)t) 

(7=1 1=1 j = l 

2 n m 

^ E EE II <t>ij Ik +^nmC v | A | i? 2 (exp(2At) - 1) 

cr=l i=l j=l 

+R(exp(Xt) - 1) , 
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which gives (38) if 

2 Ti TTb 

#>£££ Wloo=#i, (42) 

cr=l i=l j=l 

and 

1/ 9 A 2 n m \ 

4 £ T = 2A'" 1 + ^lAiff C-EEE li « I-)) • («) 

x '/II (J = l 2=1 J = l ' 



Lemma 4.4 Let T given fry wi/i A = RC V | A | (i? > R\). Then then 
there exists 4>o,R2 such that if \\ cf) ||i< 4>q, R > R2 one has 

AC 

(AT-1) 2 > — H0IU, (44) 
where C = Max(C\,C2) is given by (31) and (39). 

Proof of Lemma 4.4: Using the inequality ln(l + x) < x, x > 0, then one 
has: 

* T < ^ = — , (45) 

nmCrj I A I it nm 



which yields 

(AT - l) 2 > ( — — - ) ■ (46) 



nm — 1 x 2 



nm 



Then, if (f> and R are such that 



(nm - l) 2 C v A 
8(nm)^ (C V L A + L^J 

and 

^ 8 ,, , „ (nm) 2 „ . , , 

I A I (nm - l) 2 

Then (44) is proved. This completes the proof. ■ 

Now we have all tools for the proof of the local existence theorem. 
Proof of Theorem 3.1: Let us consider the following subset in Xt defined 
as follows: 

B T = {ip = (V^-) £ X T : Y£ > 0, II m ||x< a || (/> ||i, 

2 n m 

x - v i «>s((?i)t, y - f j sin(6»j)t) 

<T=1 1=1 J = l 

< Rexp(Xt),t G [0, T], x G 0}. 
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Let ifj l ,ip 2 e Bt, and let R> Rq = Max(Ri, R2) then from Lemma 4.3, 
there exists T such that (A^ 1 ))^ > and satisfying (38) for A = RC V | A |, 
t G [0, T]. Moreover from Lemmas 4.1 and 4.3, one has: 

II A^ 1 ) \\ Xt <\\ 4> ||i +jCa 2 || <f> \\j +Xa T \\ 4> ||i, (49) 



and 



|| A^ 1 ) - A{^ 2 ) \\ Xt < (hoa || 4> ||i +\T) I) V 1 - ^ ||x T , (50) 

where C is the constants defined in Lemma 4.4. Let 

4-C 

A O = (AT-1) 2 - T ||0||i. (51) 

Then by the previous Lemma one has Ao > 0, for || ||i< 0o- Let now ao 
be the positive quantity given by: 

(1 - AT) — -y/Ap 

ao = A 2cu\\ 1 ' (52) 

which is a solution of 

|| <p Hi +^Coo II Hi +^«oT || </> ||i= ao || ||i • 
That shows from (49) that A^ 1 ) £ Sy. Moreover 

^ Ca || ||i +AT) = 1- v/Ao < 1 

which shows from (50) that the mapping A is a contraction in a ball Bt, 
and application of the fixed point theorem gives the proof of Theorem 3.1, 
which refers to local existence. Moreover, if 4> satisfies (18), and if || <j> ||i is 
very small, then from (42) and (48), one can choose R such that (19) can 
be obtained. This completes the proof. ■ 



5 Proof of Existence for Large times 

The aim of this section consists in proving Theorem 3.2 concerning the global 
existence of solutions to the initial value problem for Eqs.(5)-(9). It will be 
proved that the solution can be extended in each interval [0,pT], for p £ N. 
As in the preceding section some preliminary results are needed: 
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Lemma 5.1 Let <fi satisfy the conditions of Theorem 3.1. Consequently, 
there exists 4> l such that if || <j> ||i< <p l , there exists a\ such the solution to 
Eqs.(5)-(9) can be extended in the interval [T, 2T] and satisfies the estimate: 

|| ip(t + T) ||i< ai || ||i, te[0,T], (53) 

p(t + T,x)<R, Vte[0,T], xeft. (54) 

Proof of Lemma 5.1: Let 

Ai = (AT-l) 2 -^Ca o ||0||i, (55) 
where ao,T, and A are fixed and depending on 4>q and Rq. Then for 

\2 



*<? =u ^m=iL) (56, 



the following Ai > holds true. 

We solve the problem (5)-(9) in [T,2T] with initial condition given by 
Y>(T;x). Then, for any t £ [0,T], one has: 

$5(t + r,x) = ^.(r,x) + jf t ^'K^(5 + r,x)exp(-A(a + r)) 

+ $£(s + T, x)|a- A^](s + T, x) exp(-A(s + T))|^)(&7.) 
Consider now the ball 

B l T = {g(t) =# + T)eI T :^>0, || ^(i + T) ||i< oi || ||i, 

2 n m 

EEE^ + T - a; - v i cos((9i)(t + T),y — Vj S m(0i)(t + T)) 

cr=l i=l j'=l 

< J Rexp(A(t + T)), tG[0,T],xG^| 
Using the same technique as in the proof of Theorem 3.1 yields: 
|| A(^)(t + T) \\ Xt < a || ||i +\ca\ \\ <j) \\j +\ ai T \\ 4> || 1} (58) 

|| A^ 1 ) - A{^ 2 ) \\ Xt < Qcai || <p ||i +\T) I) V 1 - ^ ||x T , (59) 
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when ifj l ,ip 2 e B\.. By the same arguments as in the proof of Theorem 3.1, 
and by choosing a\ given by 

the fixed point Theorem yields existence of solution in [T, 2T] . This solu- 
tion is continued in [T, 2T] and in particular, satisfies (53) and (54). This 
completes the proof of Lemma 5.1. ■ 

Proof of Theorem 3.2: The iteration process can now be applied to prove 
global existence in [0, oo[. Suppose that for small initial data, the solution 
exists and is continued [0, (p — 1)T] satisfying: 



||^(* + (r-l)T)||i<a r _i||^||i, r = l,....p-l, t€[0,T], 

and 

II ip{rT) ||i< a r _i || <j) ||i, r = 1, ....p - 1 , 
where a r are given by the following: 

°- =a (1 "2C?;||^ ' r = !,....„ -I, (61) 

A r = (AT- l) 2 -4Ca r _iW^, r = l,...p-l. (62) 

A 

and the reels Ao and ao are given respectively by (51) and (52). 

It can now be proved that the solution can be extend in [(p — l)T,pT], 
satisfying for any t € [0, T] the following: 

|| + (p - ||i<ap_i ||0 ||i, (63) 

II 4>{pT) ||i< a p _i || ||i . (64) 
Let Y> = (^) be the solution of the following problem: 

^§(t + (p - 1)T, x) = V£((p - 1)T, x) 

f t+(p-l)T 

' rf.[^,^( S ,x)exp(-A S )^(s, x) 



(p-l)T 



x |a-L[Y>](s, x)exp(-As)j^ds. (65) 
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The same arguments as in the proof of Theorem 3.1 yields 

Ca 2 

II A{^){t + (p - 1)T) || Xt < a p _ 2 || <j> ||i +— ^ || ||? +Xa p - 1 T \\ <f> ||i, 

(66) 

and 

II + (P - - ^ 2 )(t + (p - 1)T) \\ Xt 

< ^2CV_i 

II ^ Hi +Ar) II ^(t + (p - 1)T) - ^ 2 (t + (P - 1)T) \\x T ■ (67) 



The proof is complete if we choose a„_i such that 



where 



_ (l-AT)-yA^T 

" A — 2cIJTi — ' (68) 



A p _ 1 = (AT-l) 2 - 4 % 2 ll^lk, (69) 



A 

which give the solution in [(p — l)T,pT], satisfying (63)-(64). Then the proof 
is completed. ■ 



6 Simulations 

This section shows how numerical solutions can be obtained for the initial 
value problem in unbounded domains, namely within the framework of the 
qualitative analysis of Sections 4 and 5. An important technical difficulty is 
the computation of the left hand side of the equation modeling the transport 
term. The approach developed in this paper is based on the so-called split- 
ting method, where the equation is splitted into the transport part and the 
interaction term. High resolution finite differences methods for hyperbolic 
equations [17] are used. 

In this paper, we propose to solve the system of partial differential equa- 
tions by using a splitting method. As stated in [22, 31], the idea behind this 
type of approach is that the overall evolution operator is formally written 
as a sum of evolution operators for each term in the model. In other words, 
one splits the model into a set of sub-equations, where each sub-equation 
is of a type for which simpler and more practical algorithms are available. 



21 



The overall numerical method is then formed by picking an appropriate nu- 
merical scheme for each sub-equation and piecing the schemes together by 
operator splitting. 

Let us suppose that the movement occurs with a unique velocity modulus 
v. Thus, this assumption of homogeneous velocity modulus let us introduce 
the notation f?, instead of ffy After discretizing the velocity directions, 
Eq. (4) reads: 

d t ff(t,x) + v cos(^)^/f(t,x) +v sin^Ofy/f (*> x ) = J[f](t,x). (70) 
Then, we can split this equation into the following set of sub-equations: 
dtf?(t, x) + v cos(0i) 8 x f[(t, x) = 0, (71) 
d t ft (t, x) + v sin(^) d y f[(t, x) = 0, (72) 

and 

d t /f(i,x) = J[f](i,x). (73) 

Notice that Eqs. (71)- (72) are one-dimensional homogeneous transport 
equations, while Eq. (73) is just an ordinary temporal differential equation. 
We will denote by S x , S y and St the solution operators of Eqs. (71)-(73), 
respectively. For instance, S x uq is the solution of Eq. (71) for an initial 
condition uo(x). 

It is worth mentioning that one of the main advantages of this method 
is that the three new sub-equations can be properly solved by using a finite 
difference scheme. With that purpose, we discretize the time interval [0, T] 
by taking N t equally spaced points tk = kAt, k = 0, . . . , N t — 1. Now, given 
an initial condition ff = /j(0,x) for Eq. (70), we propose the following 
splitting process: 

Algorithm 
Input: fl 

For k = 0, . . . , Nt — 2 perform the following steps: 

- Solve Equation (71) with initial condition ff, obtaining a solution 

- Solve Equation (72) with initial condition S x ff, obtaining a solution 
SySxfi ■ 

- Solve Equation (73) with initial condition S y S x f^, obtaining a solution 

StSyS X ff. 
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Figure 3: Initial condition for case-study 1. A crowd localized in the bottom- 
left corner is constituted by individuals moving in five different directions, 
indicated by the arrows. The parallel arrows on the right show the exit of 
the room. 



- Set = S t S y S x fl 
Output: k = 0,...,N t -l. 

The output of the algorithm is the numerical solution of the original 
problem (70), given by f[(kAt, x) = ff, k = 0, . . . , Nf — 1, while the inter- 
ested reader is referred to [22] for a detailed analysis of the convergence of 
the method under consideration. 

In the next paragraphs we propose two case-studies in order to carry out 
some simulations, to analyze the accuracy of the model and of the proposed 
numerical scheme and to investigate the emerging behaviors. 

Case-study 1 

Let us first consider a group of individuals initially located near the corner 
of a rectangular room represented by the domain Q = [0, 1] x [0, 1] C R 2 , as 
shown in Figure 3. 

Let us suppose that all of them are moving, at time t = 0, with constant 
velocity v but with five different directions given by 9{ = in/ 8, i = 0, . . . ,4. 
Moreover, the individuals are equally distributed among these directions. In 
that moment some panic situation occurs at that corner and all of them try 
to move towards the exit, which is indicated by the parallel arrows in Figure 
3. 

Notice that those individuals who are moving in the direction 9q = are 
the most favored, as they are already moving towards the exit. However, we 
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(b) 




(c) (d) 

Figure 4: Case-study 1. Evolution of the density of the crowd for different 
times, (a) corresponds to the initial density, while (b)-(d) show its evolution. 
The parallel arrows indicate the target. 

assume that the other individuals start to change their directions according 
to a table of games considering only the sensitivity to reach the target. 
Assuming that the activity is homogeneously distributed, i.e. uq is constant, 
that A(x) = x, and taking v = 0.03 and a\ = «2 = a = 0.06 we perform 
some simulations, obtaining the results shown in Figure 4. 

We can observe that as time evolves the pedestrians are closer to the 
exit, because many of them were able to change their directions. It is worth 
mentioning that the fact of considering only the tendency to reach the target 
let the individuals correct their movements without getting confused by the 
tendency to follow the stream. 

In Figure 5 we show a top view comparing the cases a = and a = 
0.06 for a fixed time. It becomes clear that originally the five groups were 
following different targets and that the model is successful predicting the 
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(a) (b) 

Figure 5: Case-study 1. Comparison between the evolution for (a) a 
and (b) a = 0.06. 



reaction after a panic situation. 
Case-study 2 

Let us consider two different groups, labeled by Gl and G2, which are ini- 
tially moving towards two different targets inside the rectangular room, as 
shown in Figure 6. Both groups have the same velocity modulus v and we 
consider five directions of movement 0j = i = 0, . . . , 4. In this way, Gl 
and G2 are initially moving in the directions 9$ = and #4 = ir, respectively. 

Once more, we assume that they realize about the location of the exit, 
which corresponds to the direction 82 = 7r/2, and they consequently try to 
adapt their movement. Notice that in this case none of the groups originally 
have the proper direction of movement. 

The table of games takes into account that the individuals can change 
their direction with a certain probability in order to reach the exit. However, 
we consider that individuals from group Gl have better abilities, changing 
their direction with a probability a\ = 0.2, while individuals from group G2 
do it with a probability 02 = 0.1. The simulations are performed with these 
values and taking A(x) = x and v = 0.03. The results are shown in Figure 
7. 

In Figure 8 we show a top view comparing the afore-said situation with 
the case in which ct\ = a<i = 0.2. Notice about the asymmetric shape of 
the crowd in Figure 8(a), clearly given by the different abilities (e.g. visual, 
fitness) that each group exhibits. 
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X 



Figure 6: Initial condition for case-study 2. Two crowds originally moving 
in opposite directions, indicated by the arrows. The parallel arrows in the 
top show the exit of the room. 

7 Looking Ahead to Perspectives 

Various approaches to modeling crowd dynamics corresponding to different 
representation scales and related mathematical tools have been reviewed and 
critically analyzed in the survey paper [6] . The mesoscopic approach in this 
present paper has been selected, out of this overall analysis, as the most ap- 
propriate, according to the authors' bias, to capture the complexity features 
of the system under consideration. Moreover, the use of discrete velocity 
variables has been adopted to model the micro-scale state of pedestrians in 
order to capture the granular features of their dynamics [13]. 

The use of discrete velocity variables was introduced, as already men- 
tioned, to model vehicular traffic[15, 12], while the approach was further 
developed toward crowd modeling in [3]. The discretization is here applied 
to the velocity variable. In principles, it can be referred also to the space 
variable as already known for models of vehicle dynamics [16]. However, 
the technical difficulty of dealing with more than one space variables creates 
nontrivial technical problems, which do not seem yet solved at present. 

The modeling approach proposed in this paper refers to normal flow con- 
ditions, namely when pedestrians are sufficiently relaxed to behave normally. 
On the other hand, the behavioral influence of panic conditions can play an 
important role in the modification of the usual dynamics [19, 20]. This ob- 
jective can be pursued by an appropriate tuning of the model proposed in 
this paper also by taking advantage of recent studies on the behavioral dy- 



26 




Figure 7: Case-study 2. Evolution of the density of the crowds for different 
times, (a) corresponds to the initial density and the arrows on each crowd 
indicate the initial direction, while (b)-(d) show the evolution. The parallel 
arrows indicate the target. 




(a) (b) 

Figure 8: Case-study 2. Comparison between the evolution for (a) a\ = 0.2, 
«2 = 0.1 and (b) a.\ = 02 = 0.2. 
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namics of pedestrians [25, 27], and evolutive game theory [11]. Accordingly, 
the guidelines toward this important objective are given: 

• Pedestrians undergo a trend to higher values of the activity variable. 
Moreover, the dynamics over the said variable has to be necessarily 
inserted in the model; 

• The effective visibility zone becomes larger and signals from large dis- 
tance become important, while in the case of normal conditions short 
distance signals appear to be more important; 

• Pedestrians appear to be sensitive to crowd concentration, while in 
normal conditions they attempt to avoid it. Therefore different weights 
e need to be used; 

• The role of boundary conditions needs to be properly modeled ar- 
guably by means of nonlocal conditions 

Moreover an important topic to be studied is the role of control actions 
by external systems. Namely the system needs to be modeled as an open 
one. 
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